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The  results  of  our  research  in  the  first  eight  months  of  the  contract 
year  1977-78  are  summarized  in  the  three  attached  papers  which  were  pre- 
sented at  the  11th  International  Symposium  in  Rarefied  Gas  Dynamics.  Two 
of  them  will  appear  in  the  Symposium  proceedings. 

Theoretical  research  has  since  proceeded  in  two  directions: 

(1)  The  model  presented  by  Trilling  for  spherical  droplets  of 
continuous  material  is  generalized  to  other  geometrical 
configurations.  Since  ice  crystals  appear  to  grow  as 
hexagonal  platelets  it  has  appeared  interesting  to  consider 
the  growth  of  flat  circular  cylinders  (R  >>  h where  R is 
the  cylinder  radius  and  h is  its  height).  The  rate  of 
growth  on  the  sides  turns  out  to  be  sufficiently  faster 
than  on  the  circular  faces  to  explain  why  a thin  pillbox 
does  not.  grow  into  a sphere.  The  rate  of  growth  is 
fastest  near  sharp  edges  and  therefore  edges  are  not 
smoothed  out. 

(2)  Some  further  extension  of  the  statistical  mechanics  of 
nucleation  embryos  beyond  the  enclosed  paper  by  Gobbini 
appear  in  his  S.M.  Thesis.  However,  to  our  great  regret, 
Gobbini  has  now  left  this  project  to  do  his  doctoral 
research  in  another  field,  and  we  are  seeking  a replace- 
ment for  him  to  continue  the  statistical  theory. 

Experimental  research  beyond  the  results  of  the  symposium  paper  has 
necessitated  some  rather  extensive  rebuilding  of  our  experimental  apparatus. 
In  particular  we  have  been  in  the  process  of  reconstructing  the  molecular 
beam  source  chamber,  and  the  target  mount  assembly.  Extensive  experimental 
measurements  on  the  scattering  of  beams  of  water  vapor  from  ice  and  Covellite 
surfaces  have  confirmed  and  extended  earlier  results  obtained  in  our 
laboratory  by  A.  Swartzon. 
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Mass  Distribution  of  Droplets  in  Nucleation 
Leon  Trilling* 

Massachusetts  Institute  of  Technology 
ABSTRACT 


We  present  a general  model  for  the  interaction  of  a speci- 
fied distribution  of  water  droplets  with  a supersaturated 
atmosphere.  In  terms  of  an  intrinsic  nucleation  time  scale, 
the  evolution  of  the  droplet  size  distribution  is  obtained 
explicitly  as  a function  of  the  initial  state  and  of  the 
rate  of  formation  of  new  nuclei.  The  connection  between  the 
intrinsic  nucleation  time  scale  and  real  time  is  then  found 
by  solving  the  bulk  conservation  equations  for  the  atmosphere 
with  the  condensation  as  an  energy  source  term.  Examples 
are  worked  out  for  an  initially  uniform  droplet  distribution 
for  various  rates  of  nucleus  formation. 


In  a number  of  condensation  problems  of  practical  impor- 
tance, it  is  useful  to  know  the  mass  distribution  of  droplets 
or  crystallites  as  condensation  proceeds.  Depending  on  the 
order  of  that  mass,  the  problem  is  best  approached  by  means 
of  molecular  kinetics  and  statistical  mechanics  (particles  of 
up  to  a few  thousand  molecules  in  a rarefied  vapor  or  gas 
vapor  mixture)  or  as  a two  phase  continuum  non-equilibrium 
thermodynamics  problem  (solid  or  liquid  nuclei  of  order  0.1 
to  10  microns).  In  the  latter  case,  the  mean  free  path  of  the 
vapor  or  gas-vapor  mixture  may  be  large  or  small  compared  to 
the  diameter  of  the  nucleus. 

The  statistical  mechanics  of  small  crystallites  (^lO2  to 
105  atoms)  interacting  with  their  vapor  raises  two  difficult 
questions  - What  is  the  frequency  distribution  of  their 
configurations  in  phase  space?  What  is  the  energy  (particu- 
larly the  vibrational  energy)  to  be  ascribed  to  any  given 
configuration?  The  first  of  these  questions  can  now  be 
answered  fairly  well  (see  for  example  F.  Abraham  [1]  or  the 
paper  submitted  to  this  symposium  by  G.  Gobbini  [2]).  The 
answer  to  the  second  is  still  subject  to  considerable  argument, 
jand  depends  on  the  model  chosen  for  the  crystallite.  Such 
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models  require  either  a definite  set  of  rules  for  individual 
molecular  interactions  or  a convincing  way  of  relating  the 
behavior  of  small  assemblages  to  that  of  the  material  in  the 
bulk. 

If  the  nucleus  is  large  enough  (a  diameter  of  0.1  micron 
corresponds  to  some  108  atoms)  it  is  usually  represented  as 
a solid  or  liquid  sphere  whose  main  properties  (e.g.  surface 
tension,  chemical  potential)  are  those  of  the  material  in  the 
bulk.  Because  its  mass  is  large,  it  behaves  in  collisions 
with  ambient  molecules  of  vapor  or  gas  as  a body  at  rest. 

The  nucleating  particle  size  distribution  in  a non-equili- 
brium situation  is  governed  by  a supply  process,  a growth  (or 
decay)  rate  which  depends  on  the  degree  of  background  satura- 
tion, on  the  fluid  mechanics  and  surface  physics  of  inter- 
actions with  the  background,  and  on  an  exhaustion  process 
which  removes  particles  from  the  field  at  a rate  which  depends 
on  particle  size. 

We  discuss  the  nucleation  of  water  droplets  (supplied  for 
example  by  an  aircraft  jet  exhaust  or  some  other  device)  in  a 
supersaturated  atmosphere  of  specified  initial  properties, 
under  conditions  (e.g.  initial  droplet  size)  which  makes 
re-evaporation  sufficiently  rare  to  be  negligible. 


Let  the  ambient  vapor  density  be  pv  and  the  corresponding 
saturation  density  be  pv  ; then  the  saturation  ratio  is 
defined  as: 
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and  pV()  is  related  to  the  ambient  temperature  T by  the  Clausius 
Clapeyron  equation: 
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(2) 


where  L is  the  latent  heat  of  condensation  per  molecule. 


If  the  vapor  molecules  and  the  droplets  or  crystals  are 
convected  by  the  surrounding  air  motion,  then,  applying  the 
continuity  equation,  for  both  air  and  vapor,  to  eliminate  the 
divergence  of  the  velocity  vector,  we  obtain 
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i where  I is  the  total  rate  of  condensation  of  the  vapor  per 
unit  volume  and  unit  time  and  pa  is  the  air  density.  We 
write  the  energy  conservation  equation  for  the  ambient  air 
under  the  assumption  that  the  heat  released  by  condensation 
is  fully  transmitted  to  the  air  (e.g.  Pv/Pa  <<!) 

^ c„  AT  dp 


CP  il 


ST  + LI 


This,  with  the  (ideal  gas)  equation  of  state  for  air  eliminates 
the  temperature  and  density  from  (1-3);  finally  eliminating 
dp«0/dt  by  use  of  (2),  we  obtain: 


1 dS 
1-S  dt 


d(logpa) 

dt 


(1  + ^ (W>2) 


This  is  the  bulk  equation  which  connects  the  saturation  ratio 
I history  to  the  dynamics  of  the  bulk  air  motion  and  to  the 
(total  vapor  condensation  rate.  A similar  derivation  may  be 
found  in  [ 3 ] . 

If  we  now  represent  the  condensed  phase  as  a distribution 
of  spherical  clusters  f(r,t)  where  f(r,t)  is  the  number  of 
(clusters  per  unit  volume  whose  radius  is  in  the  range 
(r,r  + dr),  then  the  condensation  rate  I is 


r2f (r,t)  dr 

dt 


The  continuity  equation  for  the  cluster  size  density  in 
(r,t)  space  is 


h «£  > 


where  ^(r.t)  is  the  rate  per  unit  volume  per  unit  time  at  which 
new  nuclei  of  mass  (r,  r+dr)  appear  in  the  fluid. 

i 

Finally,  we  balance  the  energy  transferred  to  a droplet 
by  collision  with  ambient  gas  molecules  (or  by  convection) 

.and  the  energy  exchanged  in  the  condensation  process,  and 
obtain  a relation  of  the  form: 


fr  “ F(r,s)  - sv(r)  (8) 

dt  | 

When  the  time  scale  of  collisions  is  short  compared  to  the  time 
scale  for  droplet  growth  as,  in  the  Stokes  flow  regime,  Fukuta 
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• and  Walter  4,  and  independently  Sedunov3 , obtain 
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Here  pa  is  the  viscosity  of  the  air 

Dav  is  the  diffusion  coefficient  of  water  vapor  in  air 

mv,ma  are  the  molecular  masses  of  vapor  and  air 
respectively 

h is  a characteristic  length,  essentially  a Knudsen  layer 
thickness 


ycp/k 
V/PaDav 


Pr  is  the  Prandtl  number  of  air  Pr 
S is  the  Schmidt  number  of  air  S 
o is  the  energy  accommodation  coefficient  for  air  on  water 
3 is  the  condensation  coefficient 

The  equivalent  result  for  the  free  molecule  regime  is: 
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Note  that  i.a,  ig  and  h are  lengths  of  the  order  of  the  mean 
|free  path  of  the  gas;  A,  B are  numbers  of  order  unity.  We 
expect  (9)  to  hold  since  the  Reynolds  number 


Pa 


r dr/dt 


is  of  order  paS/PL  which  is  well  below  unity  when  the  radius 
of  the  droplet  is  larger  than  the  mean  free  path.  In  air  at 
sea  level,  the  mean  free  path  of  air  is  of  order  0.1  microns; 
one  would  therefore  use  (9)  for  all  droplets  considered 
there.  But  at  an  altitude  of  50  km,  the  mean  free  path  is 
100  microns  and  formula  (10)  is  more  appropriate. 
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Substituting  (8)  into  (9)  or  (10)  and  introducing  the 
new  coordinates 
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j we  obtain 
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V(r,t)  = <Ji  (o ,x) 


(12) 


which  is  integrated  formally  to  give 

x 

g(a,x)  = g0(o-t)  + { 


iKC+a-x,  O d£ 


(13) 


i where  jLflXg).  is  the  initial  distribution  of  droplet  sizes, 
vo  (a) 

Thus,  in  the  absence  of  a source  of  nuclei,  g shifts  to 
! the  right;  the  distribution  function  keeps  its  shape  but  is 
applied  to  larger  droplets;  the  droplets  in  the  system  grow 
■and  the  distribution  function  scale  is  modulated  by  the 
factor  v„/v  (see  11c)  which  is  generally  greater  than  unity. 
The  time  scale  t is  not  real  time;  it  is  related  to  the  satur- 
ation ratio  by  (11a)  and  as  S->-0,  t approaches  a finite  value 
while  t-*00;  when  the  droplets  have  condensed  all  the  available 
vapor,  they  grow  no  further. 

If  a source  term  is  specified  [e.g.  V (r,t)],  the  form  of 
ji|i(o,T)  is  not  known  since  (11a)  cannot  be  inverted  without 
(solving  (5);  but  its  role  is  clear;  it  provides  new  nuclei 
iwhich  compete  with  the  original  ones  in  g0(o)  for  available 
ivapor  and  therefore  inhibit  the  growth  of  large  droplets  by 
(drawing  some  vapor  to  new  nuclei.  As  a result,  the  condensa- 
tion process  is  completed  sooner  in  real  time. 


To  complete  the  formal  solution  of  the  problem,  we  substi- 
tute (13)  into  (5)  transformed  to  the  t coordinate: 
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j Remembering  the  definition  of  S,  we  eliminate  pv  from  (14), 


since  1-S  = 
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Let  Ij  denote  the  integral  of  the  homogeneous  part  of  the 
solution  and  I2  that  of  the  non-homogeneous  part.  In  the 
special  case  of  an  atmosphere  at  rest  d log  pa/dt  vanishes, 
and  in  terms  of  Ij  and  I2, 
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Integrating  (14b) 
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where  is  the  initial  supersaturation.  Whenever  the  expo- 
■ nent  of  the  exponential  is  positive,  S < S^,  and  as  the  expo- 
nent increases,  S decreases  monotonically  toward  zero. 

The  real  time  variable  t is  obtained  by  integrating  (11a) 
when  S(x)  is  known.  In  particular,  saturation  occurs  when  S 
vanishes,  or 


[Ij(t)  + i2(r)]dx= 
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since  .‘l  «1  for  most  problems  of  practical  interest,  (16) 
Y Pa 

may  be  approximated  as 
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where  pv^  is  the  initial  ambient  vapor  pressure  and  Pv0 

saturation  vapor  pressure  at  the  (constant)  ambient  tempera- 
ture. The  integrals  Ii(x)  and  I2(x)  represent  the  contribu- 
tions of  the  growth  of  original  droplets  and  added  source 
droplets  respectively  to  the  depletion  of  the  available  vapor 
isupply.  The  integral  Ij  can  be  computed  explicitly  when  the 
original  distribution  of  droplet  size  f0(r)  and  the  growth 
law  (8)  are  known.  The  integral  I2  includes  < p(a,x)  and  cannot 
ibe  computed  directly  from  ip(r,t).  But  a step  by  step  numerical 
integration  where  12  is  evaluated  approximately  over  the  nth 
step  by  using  the  value  of  Sn_^,  gives  a good  approximate 
result  in  a straightforward  manner. 
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For  Stokes  flow  (9),  a convenient  non-dimensional  time  x 
based  on  the  rate  of  growth  of  a reference  droplet  of  radius 
r0  is  introduced: 


pLr0(Ar0  + Bh) 


(18a) 


For  a 1 micron  droplet  at  sea  level,  one  unit  in  x corresponds 
, to  100  seconds.  Similarly  a scale  for  the  distribution  func- 
I tion  is  introduced: 


-47TPLrofo 


(18b) 


That  is  essentially  the  ratio  of  mass  of  condensed  material 
to  the  mass  of  vapor  at  saturation  conditions  in  the  same 
volume.  For  a set  of  droplets  of  size  r0,  the  integral 
j/l1dt  is  evaluated  as 


Ij  dt  = n(T  + x2) 


(19a) 


land  the  integral  / I2dx  becomes: 
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jwhere  £ is  the  non-dimensional  source  intensity  defined  in  a 
■manner  comparable  with  r)  and  t is  a non-dimensional  "real 
:time"  compatible  with  X.  The  effect  of  the  source  term  is 
i shown  on  Figures  1,  2 by  comparing  the  real  time  histories 
lof  the  condensation  for  n = 1,  5 = 0,1,5  for  an  initially 
uniform  set  of  droplets  r0. 
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ABSTRACT 

Homogeneous  nucleation  depends  on  the  free  energy  of  con- 
densational embryos,  an  important  part  of  which  depends,  in 
its  turn,  on  the  embryos'  natural  frequencies  of  vibration. 
Mathematical  models  previously  proposed  for  these  vibrations 
are  discussed,  then  a new  one  is  proposed,  which  envisages 
both  compressional  and  distortional  vibrations.  This  model 
is  applied  to  small  cubic  crystallites,  and  the  resulting 
values  of  free  energy,  obtained  by  means  of  numerical  compu- 
tation, are  presented  and  discussed. 

It  is  well  known  that  the  liquefaction  or  sublimation  of  a 
gas  is  severely  delayed  when  there  are  no  pre-existing  nuclei 
on  whose  surface  condensation  may  begin.  Under  such  condition^ 
the  supersaturated  gas  is  believed  to  build  embryonic  aggre- 
gations, or  clusters.  For  smaller  clusters,1  free  energy 
increases  with  size,  and  so  these  clusters  tend  to  decay. 

Above  a critical  size,  however,  free  energy  decreases  with 
size.  When  supercritical  clusters  are  formed,  they  tend  to 
grow,  and  condensation  finally  occurs. 

The  free  energy  of  a cluster  made  of  n molecules  is: 

fn  = -kT  ln(Qne  mn/Nn) 

where  m:  total  number  of  molecules  in  the  enclosure 
Nn:  number  of  clusters  made  of  n molecules 

and:  Qn  = (n!)-1  I exp  (-Hjj/kT)  dx  d"p 

where  Hn  is  the  Hamiltonian  of  the  cluster,  x and  "p  its 
Lagrangian  position  and  momentum  coordinates.  If  we  suppose 
the  clusters  to  be  crystallites,  with  their  molecules  oscil- 
lating around  invariable  lattice  positions,  Hn»  and  therefore 
Qn  and  fn,  can  be  decomposed  into  the  sum  of  terms  due  to 
translational,  rotational,  vibrational  and  binding  energy. 

The  vibrational  contribution  to  free  energy  can  be  expressed, 
if  the  temperature  is  not  too  low,  by  the  formula: 

h 3 n-6 

(1)  f_  = (3n-6)kT  In  tf  + kT  In  n w 

n»  V1D  i=i  nl 
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where  the  u)n^'s  are  the  3n-6  natural  frequencies  (in  rad/sec) 
of  a crystallite  of  n molecules.  It  has  been  pointed  out2 
that  the  "replacement"  term  of  free  energy  lies  within  this 
component . 

An  expression  for  the  intermolecular  potential  is  neces- 
sary in  order  to  find  the  natural  frequencies.  Two  different 
kinds  of  expressions  have  been  proposed  so  far,  to  my  know- 
ledge. One  is  a potential  depending  only  on  the  distance 
between  nearest  neighbors.  If  the  problem  is  linearized,  this 
means  modeling  the  crystallite  as  mass  points  coliegated  by 
ideal  coil  springs.  A serious  objection  must  be  raised  against 
this  model:  it  works  for  a f.c.c.  lattice  - to  which  it  was 

q 

applied  by  Hoover  et  al.  - but  would  not  work  for  most  other 
lattices.  In  most  crystalline  lattices  non-rigid  deformations 
are  possible,  which  leave  the  distance  between  nearest  neigh- 
bors unchanged;  the  model  provides  no  restoring  force  for 
such  deformations. 

The  directionality  of  bonds,  which  holds  these  lattices 
together,  was  taken  into  account  in  the  potential  proposed  by 
Stillinger  and  Rahman  for  liquid  water4’  . The  water  molecule 
was  modeled  as  a rigid  tetrahedron  with  four  electrical 
charges  at  the  corners,  and  the  potential  was  the  sum  of  a 
term  depending  on  the  distance  between  centers  and  of  the 
electrostatic  terms.  The  motions  allowed  were  compressions 
and  librations  - the  latter  being  rigid  rotations  of  a mole- 
cule with  respect  to  its  neighbors.  If  a model  like  this 
were  applied  to  a crystallite,  it  would  account  for  the 
stability  of  any  kind  of  lattice,  but  would  probably  envisage 
unrealistic  vibrational  modes.  The  assumption  of  molecular 
rigidity  is  too  strict.  In  the  very  case  considered  in 
ref.4,5#  water,  the  equilibrium  angle  between  valence  orbitals 
changes  of  about  5%  from  the  water  vapor  molecule  to  the 
perfect  ice  crystal.  It  is  only  reasonable  to  assume  that 
oscillations  of  this  angle  around  its  equilibrium  value  occur 
in  crystallites.  On  the  other  hand,  allowing  librational 
freedom  to  molecules  is  probably  too  liberal  an  assumption  for 
crystalline  aggregates,  where  bonds  are  complete  and  fairly 
tight,  while  it  seems  correct  for  a liquid,  where  bonds  are 
incomplete,  strained  and  easily  broken. 

For  these  considerations,  a different  linear  model  of 
crystalline  lattices  is  proposed  here.  It  is  supposed  that 
molecules  are  not  rigid,  but  that  their  valence  orbitals  always 
stay  oriented  toward  the  centers  of  nearest  neighbors.  Libra- 
tions, i.e.,  rigid  rotations  of  a molecule,  are  therefore 
excluded.  The  displacement  of  a molecule  from  its  equilibrium 
position  will  cause  a change  both  of  the  distance  between 
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centers  and  of  the  angle  between  orbitals.  Both  are  opposed 
by  a restoring  force,  so  that  compressional  and  distortional 
vibrations  result.  Let  us  call  Kc  and  Kd  the  elastic  constants 
of  compression  and  distortion,  respectively;  dij  the  distance 
between  the  centers  of  the  i-th  and  the  j-th  molecule;  aijk 
the  angle  with  vertex  in  the  center  of  the  i-th  molecule  and 
sides  going  through  the  centers  of  the  j-th  one  and  of  the  k-th 
one;  d°*  = d and  a°jLjk  the  respective  equilibrium  values. 

The  small-displacement  potential  will  be: 


2 i1^  Kc(dij"dij)2  + 2 ji  ^ j,k  Kd(aijk"aijk) 


2d2 


where  the  symbol  ~i  has  been  used  to  mean  a summation  over  the 
nearest  neighbors  of  the  i-th  molecule,  and  d has  been  intro- 
duced to  give  ^ the  same  dimensions  as  Kc. 

This  potential  has  been  used  to  compute  numerically  the 
vibrational  component  of  the  free  energy  of  cubic  crystallites 
with  a simple  cubic  lattice.  Crystallites  of  8,  27,  64  and 
125  molecules  have  been  considered.  The  27  molecules  cube  is 
shown  in  fig.  1.  The  force-constant  matrix  F is  found  from 
the  potential  <J>: 


Fik  ■ 1/®  32<t/3xj3xk 


m being  the  mass  of  each  molecule  and  xi,  xfc  two  of  the  3n 
coordinates.  The  eigenvalues  of  F are  the  squares  of  the  m^'s 
(see,  e.g.,  ref.1  , p.110).  The  product  of  the  u>|'s  (including 
the  six  zero  values  competing  to  rigid  translations  and  rota- 
tions of  the  whole  crystallite)  is  det  F.  Obviously,  then, 
detF  = 0.  With  a simple  mathematical  transformation  a matrix 
F'  can  be  found,  which  has  all  the  eigenvalues  of  F except 
the  zeroes.  The  determinant  of  F',  equal  to 
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has  been  computed  numerically  for  several  values  of  the  ratio 

q * / K^ / Kc  , in  order  to  express  the  results  in  a convenient 
way,  let  us  rewrite  eq.  (1)  as  follows: 


(2) 


fn,vib  = [~(3n-6)ln|+  K]  kT 


with:  0 
K 


Uc  h/k  as  a characteristic  temperature  and: 


In 


3r6 

i-=i 


[u>ni/flc]  • 


The  computed  values  of  K are  shown 
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In  table  1. 
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K 


q 

n * 8 

n m 

n = 64 

n - 125 

.05 

-7.8944 

-25.752 

-45.658 

-66.166 

.10 

-3.6467 

-10.463 

-13.500 

- 9.784 

.15 

-1.0708 

- 0.803 

8.094 

30.439 

.20 

0.8475 

6.716 

25.784 

64.706 

.25 

2.4209 

13.127 

41.420 

95.667 

.30 

3.7845 

18.858 

55.718 

124.325 

.35 

5.0073 

24.008 

69.019 

151.159 

.40 

6.1278 

29.009 

81.508 

176.448 

.45 

7.7697 

33.610 

93.304 

200.378 

.50 

8.1479 

37.962 

104.489 

223.089 

TABLE  1 


Computed  values  of  factor  K of  eq.  (2)  for 
four  cubic  crystallites  (8  to  125  molecules) 
and  several  values  of  q = /Kj/K,, 


It  is  customary  (as  in  ref.3)  to  interpret  these  results 
in  terms  of  excess  entropy.  An  Einstein  frequency  we  is 
defined  as  the  frequency  with  which  an  internal  molecule  would 
oscillate  if  the  other  molecules  were  fixed.  In  our  case 

we  - /2Q3  + 16  Qj  , where  * Ac/m,  = Aj/m  . The  excess 
entropy  per  molecule  is  then  defined  as: 


(3) 


- In  n 
n i=l 


an-6 


w. 


Notice  that,  if  Te  = hwe/k  is  defined  as  the  Einstein  tempera- 
ture, the  vibrational  contribution  to  free  energy  is: 

fn,vib  " ~kT  t(3n~6)  ln  T/Te  + nSe]. 

Figure  3 shows  the  computed  values  of  Se  versus  n.  For  each 
given  value  of  q,  Se  behaves  as  expected;  it  decreases  with  n, 
as  surface  molecules  become  less  and  less  important  than  bulk 
molecules;  the  downturn  for  small  values  of  n is  due  to  the 
"replacement"  entropy  - that  is,  substantially,  to  the  fact 
that  six  frequencies  are  missing.  What  is  remarkable  is  the 
strong  dependence  on  the  ratio  q.  For  q ■ 0,  Se  goes  to 


i 


57 

5 


infinity.  Geometrical  considerations  show  that,  in  the  n- 
molecules  cube,  there  are  3(n2/3-  2)  independent  non-rigid 
deformations  that  involve  only  distortion  of  angles,  without 
compression  of  distances  between  nearest  neighbors.  Conse- 
quently, 3(n2/3  - 2)  natural  frequencies  are  a function  of 
Kd  alone  and  go  to  zero  when  q does,  whereas  we  does  not.  To 
avoid  this  inconvenience,  we  can  introduce  another  Einstein 
frequency,  depending  on  alone:  = 4ft,j  and  then  define  a 

• modified  excess  entropy: 

. In  [oj3(n-n2/3)  w.3(n2/3  -2)  ] 

w e " n 

In  II  ojj 

i=i 


which  takes  into  account  the  existence  of  purely  distortional 
and  mixed  (compressional-distortional)  natural  modes. 

Figure  2 shows  Sg  versus  n.  The  dependence  on  q remains, 
but,  at  least,  the  singularity  is  avoided.  The  results  are 
regrettably  too  few  to  show  a clear  law  in  the  dependence  of 
Sg  on  q and  n,  but  the  general  behavior  is  correct. 

It  is  apparent  that  the  results  of  a model  involving  two 
elastic  constants  are  more  difficult  to  interpret  than  those 
of  the  model  of  ref.3  (see  their  interpretation  in  ref.1, 
p.130  ff.).  Yet,  a two-constant  model  is  necessary  unless  we 
restrict  our  investigations  to  the  narrow  case  of  f.c.c. 
lattice.  Hopefully,  the  model  proposed  here  will  prove 
instrumental  in  broadening  the  field  of  data  on  the  free  energy 
of  crystallites. 
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Time-of-Fl ight  Analysis  by  the  Pseudorandom  Chopping  Technique: 

Calibration  and  Application  to  Mass  Spectrometry 

+ ++ 

David  D.  Dreyfuss*,  Robert  B.  Doak,  Harold  Y.  Wachman 

A time-of-f light  apparatus  utilizing  a cross  correlation  chopper 
for  beam  modulation  has  been  calibrated  against  a Maxwellian  stream, 
consistent  and  reliable  scheme  for  determining  the  time  origin,  to 
i p sec  for  any  arbitrary  time-of-f 1 ight  wave-form  has  been  devised. 
It  is  also  shown  that  the  need  for  deconvolution  over  shutter  function 
is  virtually  eliminated,  since  with  this  scheme  it  is  possible  to  work 
in  principle  with  an  arbitrarily  narrow  width  of  unit  slot.  A 
relatively  simple  and  straightforward  analytical  procedure  for  decon- 
voluting  the  modulated  signal  is  described. 
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1.  Introduction 

In  1966  we  reported  [1,  2]  from  this  laboratory,  on  the  calibration 
and  performance  of  a time-of -flight  (TOF ) apparatus  (fig.  1)  in  which  the 
shutter  consisted  of  a single  rotating  disk  having  four  narrow  slots  of 
equal  widths  located  at  90°  intervals  around  the  rim.  The  detector  con- 
sisted of  an  electron  beam  ionizer  which  (ideally)  produced  a sheet  of 
electrons,  which  ionized  the  molecular  stream  at  a plane.  The  ions  were 
drawn  out  of  the  stream  and  collected  on  an  electron  multiplier. 

There  are  several  limitations  to  this  scheme.  The  time  resolution  of 
the  apparatus  is  determined  by  the  ratio  of  the  open  to  closed  arc  segments 
around  the  disk  (the  chopping  period,  which  is  the  sum  of  these  two  lengths, 
is  fixed  by  the  need  to  prevent  overlap  of  fast  molecules  from  one  pulse 
and  slow  molecules  of  the  previous  pulse).  To  obtain  good  time  resolution, 
one  uses  a narrow  slot.  It  is,  in  principle,  possible  to  deconvolute  from 
the  effects  of  a wide  slot  (shutter  function),  but  the  mathematical  manip- 
ulations involved  are  very  sensitive  to  the  magnitude  of  the  noise  in  the 
data.  In  practice,  if  the  shutter  function  is  broad  enough  to  distort  the 
TOF  waveform,  then  the  best  that  deconvolution  provides  are  the  lowest 
order  moments  of  the  velocity  distribution  represented.  Usually,  in  the 
past,  it  was  found  necessary  to  sacrifice  information  on  velocity,  and  use 
a wide  enough  slit  to  transmit  a sensible  signal  in  each  pulse. 

The  ionization  detector  has  its  own  limitations.  An  electron  beam 
ionizer  is  a universal  detector,  hence  is  usable  for  all  molecular  species 
in  the  test  stream.  As  a consequence,  however,  background  molecules  are 
detected  as  well  (in  particular,  with  permanent  gases,  all  those  stream 


molecules  which  strike  closed  portions  of  the  disk  are  detected  os  back- 
ground). Statistical  fluctuations  in  this  background  are  usually  the 
major  source  of  noise  in  the  TOP  data.  Other  detector  problems  include 
the  finite  extent  of  the  ionization  region  along  the  beam  direction  and 
non-zero  ion  travel  time  between  the  instant  of  ionization  and  detection. 
The  extent  of  the  ionization  region  does  not  present  a serious  problem 
because  it  has  only  a second  order  effect  on  the  TOF  distribution  (the 
finite  shutter  function  is  a first  order  effect).  We  have  shown  [3] 
that  a considerable  flight  path  through  the  ionization  region  (up  to 
about  252  of  the  total  flight  path)  can  be  used  before  appreciable  dis- 
tortion in  waveform  occurs,  and  we  have  taken  advantage  of  this  to 
improve  the  efficiency  of  ionization.  As  to  the  effects  of  the  ion 
travel  time,  it  can  be  made  insignificant  with  sufficient  draw-out 
potentials  or  alternatively  it  can  be  accounted  for  in  interpreting  the 
data. 


2.  Pseudorandom  Chopping 

An  alternative  approach  to  molecular  beam  chopping,  which  avoids 
many  of  the  problems  cited  above  is  to  use  a pseudorandom  chopping 
scheme  (fig.  2).  This  technique  was  originally  developed  in  connection 
with  thermal  neutron  studies  [4-7],  and  has  only  recently  been  used  with 
neutral  molecular  beams  [8,  9].  In  comparison  to  the  scheme  described 
above  a larger  number  of  pulses  of  molecules  are  produced  over  each 
revolution  of  the  disk.  As  a consequence  several  pulses  are  produced 
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within  each  chopping  period,  which  overlap  at  the  detector.  By  spacing 
the  pulses  "randomly"  within  the  period,  it  is  possible  by  a simple 
arithmetic  procedure  to  undo  the  overlap  produced  [3,  9],  It  is  possible, 
in  fact,  to  attain  a duty  cycle  of  50%  for  such  a chopper  independent  of 
individual  slot  width.  With  this  chopper,  the  time  resolution  is  determined 
by  the  length  of  the  pseudorandom  chopping  sequence,  where  "length"  counts 
the  total  number  of  unit  slot  widths  (note  that  several  adjacent  locations 
may  be  open,  thus  forming  a single  larger  slot  which  is  treated  as  a set 
of  unit  slots).  The  deconvolution  procedure  recovers  the  waveform  which 
would  be  obtained  from  a chopper  with  a single  unit  slot  per  period.  While 
it  is  conceivable,  in  principle,  to  improve  time  resolution  arbitrarily  by 
taking  longer  and  longer  sequences,  in  practice,  limitations  in  electronic 
response  time  and  the  need  to  have  unit  slot  width  compatible  with 
molecular  beam  width  precludes  this  extension.  Avoiding  these  problems  by 
using  larger  diameter  disks  will  cause  other  experimental  difficulties. 
(Given  an  electronic  chopping  scheme  of  some  sort  (e.g.  chopping  laser 
excitation  of  a beam),  this  last  limitation  may  not  apply,  since  chopping 
could  be  entirely  independent  of  collimation,  and  even  better  time 
resolution  might  be  obtainable.) 

3.  Calibration 

Although  the  pseudorandom  chopping  technique  has  been  demonstrated 
for  molecular  beams,  to  cur  knowledge,  no  calibration  of  a system  using 
the  technique  has  been  reported.  We  have  performed  such  a calibration 
against  inert  gas  Maxwellian  beams.  A schematic  diagram  of  the  apparatus 
is  shown  in  figure  3.  A room  temperature  gas  reservoir  at  a pressure  on 
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the  order  of  10  Pa  supplies  a beam  through  a small  hole  in  a thin  stainless 
steel  wall.  In  operation,  conditions  were  such  that  adequate  signals  were 
achieved  for  a source  Knudsen  number,  Kn  ^ 1 . Flow  at  this  value  of  Kn  is 
perhaps  not  fully  molecular  (one  would  like  to  use  Kn  = 100),  but  it  appears 
adequate  for  velocity  distribution  measurements  from  comparison  of  velocity 
distributions  for  a series  of  decreasing  source  pressures  which  indicates 
no  measurable  change  in  distribution  when  values  of  Kn  £ .1  are  used. 

Some  examples  of  calibration  runs  are  shown  in  figure  4.  To  compare 
the  data  with  the  known  velocity  distribution  of  the  gas  in  the  source, 
a least  square  fit  of  a Maxwellian  velocity  distribution  to  the  data  was 
made.  Source  parameters  such  as  molecular  weight  and  temperature  which  are 
known  or  independently  measurable  were  held  fixed.  Fitting  parameters  were 
a baseline  position  (baseline  information  is  lost  experimentally  because  the 
signal  is  superposed  on  a large  DC  background),  and  an  overall  amplitude 
factor.  The  distribution  is  linear  in  these  quantities,  so  the  fitting  is 
straight-forward  and  well -defined.  It  avoids  the  hazards  of  peak  height 
matching  in  that  it  weights  all  data  points  equally,  and  can  be  used  for 
comparison  with  moderately  noisy  data.  Also  it  is  readily  extended  to  gas 
mixtures,  which  should  have  a velocity  distribution  consisting  of  a 
weighted  sum  of  two  different  Maxwellians. 

4.  Time  Origin  Considerations 

One  parameter  important  to  comparison  of  the  data  with  Maxwellian 
theory,  which  has  not  so  far  been  mentioned,  is  the  location  of  the  time 
origin  in  the  TOF  waveform.  Physically,  this  is  the  time  at  which  all 


molecules  represented  by  the  TOF  waveform  passed  the  chopper.  Cecause  of 
the  sharp  leading  edge  of  the  Maxwellian  TOF  waveform,  small  differences  in 
time  origin  (less  than  one  channel  (slot)  width  or  15  psec)  are  readily 
detected.  In  making  a calibration  it  is  especially  important  that  the  time 
origin  be  fixed  independently,  and  not  determined  by  matching  peak  locations 
or  the  leading  edge  or  by  some  similar  scheme. 

In  principle,  the  time  origin  is  easy  to  determine.  Using  a photo- 
trigger reference,  whose  position  relative  to  the  beam  position  is  measured, 
time  origin  can  be  determined  in  terms  of  some  constant  number  of  channels. 

Time  delays  associated  with  either  ion  drawout  time  or  electron  signal 
processing,  are  essentially  constant,  and  may  also  be  assessed.  It  is 
difficult,  however,  to  measure  these  delays  independently  with  sufficient 
precision.  It  seems  best  to  measure  them  for  the  composite  apparatus.  Of 
the  several  schemes  attempted,  the  following  appeared  the  most  satisfactory. 

Using  the  complete  system  with  molecular  beam  and  all  processing 
electronics  in  place,  electronic  delays  were  determined  from  background  gas 
detection  without  the  chopper.  For  this  process  the  electron  accelerating 
grid  in  the  ionizer  was  pulsed  electronically,  and  the  time  delay  from  the 
electronic  pulse  to  signal  rise  was  measured.  To  determine  geometric 
timing  relationships,  a series  of  TOF  spectra  were  measured  for  a constant 
intense  signal.  (The  stream  was  not  necessarily  Maxwellian:  we  used  a 
fairly  high  pressure  f^O  beam  for  which  we  could  obtain  an  excellent  signal- 
to-noise  ratio.)  These  TOF  waveforms  were  produced  by  rotating  the  chopper 
both  "forward"  and  "backwards"  (the  relative  location  of  the  photo-trigger 
and  molecular  beam  define  directions).  Then  the  relative  location  of  some 
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known  feature  (e.g.  peak  position)  for  the  two  resulting  waveforms  were 
compared  to  give  a measure  of  the  geometric  timing  relationship  which  is 
independent  of  electronic  delays.  Improved  precision  in  the  measurement 
can  be  obtained  by  making  the  measurement  at  several  chopping  frequencies. 
We  have  been  able  to  locate  the  time  origin  to  within  about  ±1/3  channel 
(typically  ± a few  microseconds  for  our  operating  chopping  period  of 
1-2  msec  and  pseudorandom  sequence  of  length  103)  using  these  techniques. 

5.  Applications  to  Hass  Spectrometry 

One  possible  application  of  the  improved  resolution  TOF  system  we 
have  developed  is  to  mass  spectrometry  (see  fig.  5).  As  mentioned  above 
in  connection  with  the  fitting  of  theoretical  curves  to  the  data,  it  is 
possible  to  obtain  a measure  of  the  relative  amounts  of  gases  in  a gas 
mixture  by  comparing  amplitude  of  the  fitting  Maxwellians.  The  advantage 
of  the  technique  is  that  it  is  non-destructive,  in  the  sense  that  ioniza- 
tion takes  place  only  after  the  necessary  information  (time  of  flight 
between  two  fixed  points)  has  been  obtained.  There  are  difficulties, 
however,  because  the  velocity  distributions  for  the  species  to  be  separated 
must  be  known  independently.  This  information  is  available,  a priori  only 
for  free  molecular  flows.  If  Maxwellian  velocity  distributions  must  be 
used,  only  a few  fairly  widely  spaced  molecular  v/eights  can  reasonably  be 
separated  since  the  separate  TOF  peaks  are  broad  and  overlapping.  In  some 
experimental  situations,  for  example,  detection  of  dimers  and  higher  order 
clusters  in  near  condensing  systems,  the  advantage  of  non-destructive 
sensing  may  outweigh  the  disadvantages  cited. 


REFERENCES 


"5 


1.  P.  B.  Scott,  "Molecular  Beam  Velocity  Distribution  Measurements" 

M.I.T.  Fluid  Dynamics  Research  Lab.  Report  65-1 , Cambridge,  MA,  ( 1 965 ) . 

2.  P.  B.  Scott  et  al.,  in  RGD  (C.L.  Brundin  ed.)  p.  1353,  Academic  Press 
New  York,  (1967). 

3.  David  D.  Dreyfuss,  "Neutral  Time  of  Flight  Mass  Spectrometry," 

SM  Thesis,  M.I.T. , Cambridge,  MA  (1977). 

4.  T.E.  Stein,  et  al.,  J.  Nuc.  Energy  A/B  16,  499  (1962). 

5.  K.  Skold,  Nucl.  Instr.  and  Meth.  63^  114  (1968). 

6.  A.  Virjo,  Nucl.  Instr.  and  Meth.  .73,  189  (1969). 

7.  W.  Glaser  and  F.  Gompf,  Nukleonik  1_2,  153  (1969). 

8.  V.  L.  Hirschy  and  J.  P.  Aldridge,  RSI  42,  381  (1971). 

9.  H.  D.  Meyer,  Max  Planck  Inst.,  Gottingen  Bericht  113  (1974). 


0 


FIGURE  3 
The  Apparatus 


FIGURE  5 


Sample  TOF  Waveform  for  a H^/Ar  Gas  Mixture.  Solid  Curve  is 
a least-squares  fit  of  a sum  ot  two  Maxwell ians  to  the  data  (fitting 
parameters:  baseline,  amplitudes). 


